Signal integration enhances the dynamic range in neuronal systems 
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The dynamic range measures the capacity of a system to discriminate the intensity of an external 
stimulus. Such an ability is fundamental for living beings to survive: to leverage resources and 
to avoid danger. Consequently, the larger is the dynamic range, the greater is the probability of 
survival. We investigate how the integration of different input signals affects the dynamic range, 
and in general the collective behavior of a network of excitable units. By means of numerical simu- 
lations and a mean-field approach, we explore the nonequilibrium phase transition in the presence 
of integration. We show that the firing rate in random and scale-free networks undergoes a discon- 
tinuous phase transition depending on both the integration time and the density of integrator units. 
Moreover, in the presence of external stimuli, we find that a system of excitable integrator units 
operating in a bistable regime largely enhances its dynamic range. 



I. INTRODUCTION 

A system operating in the vicinity of a critical state 
can present several advantages. For instance, hair cells 
of the auditory system poise themselves close to a Hopf 
bifurcation and in neuronal systems it has been pro- 
posed to provide optimal solutions for sensory stimuli 
detection (si the transmission and storage of infor- 
mation [3,1a], and computational capabilities @- These 
results motivated discussions of how the brain can, if it 
does, operate in a critical state and whether it could be 
due to self-organization arguments Q or by evolution- 
ary reasons [8(. Neural systems operating in a critical 
state also provide an alternative explanation of how the 
brain integrates the activity of distant regions [3]. In 
the critical regime, the correlation length diverges and 
neurons from different areas can effectively share infor- 
mation. Based on these arguments and on experimental 
evidences @ , it has been suggested that the brain should 
be tuned around a critical point of a second-order phase 
transition to efficiently process information 0, ■ 

Excitable media have been proved to serve as excellent 
stimulus intensity processors. Their fundamental nonlin- 
ear interactions of excitable waves confer a great capacity 
to compress several decades of stimulus intensity inputs 
into a single decade of firing rate output [llj |. This ca- 
pability, which has also been proposed to be the main 
function of neuronal active dendrites [HI, is robust for 
different networks 0, [l2l - [l4| . In many contexts, such 
as gene regulatory networks [l5j], and neuronal [l6j] and 
social systems [17] . the typical elementary unit dynam- 
ics results from the integration of neighbor contributions. 
In neuroscience, it remains a fundamental open problem 
to understand how a singular membrane potential out- 
put is generated by the convergence of complex spatio- 
temporal synaptic integration [TJ [l2|, [l8[ . To accrue for 
this difficulty, neurons present a myriad of active chan- 
nels [l9| , dendritic structures (even within the same neu- 



ron type (20[), and temporal integration modes. For ex- 
ample, the efficacy of the presynaptic neurons is largely 
variable, and neurons might require up to hundreds of 
excitatory postsynaptic potentials to spike fl6j . 

In this letter we demonstrate that integration of ex- 
citable units is a central element to shape the dynamics 
of the system: The nonequilibrium phase transition, be- 
tween the resting and the self-sustained configurations, 
switches from a continuous second-order to a discontin- 
uous first-order transition. Along with this discontinu- 
ity, a history-dependent bistable phase emerges. In this 
phase, the input-output response changes and the dy- 
namic range is strikingly enhanced. We show the gen- 
erality of the result with respect to the network topol- 
ogy, the integration time window, and the number of in- 
put signals needed to fire. Moreover, we point out how 
the presence of a bistable phase changes the paradigm 
of maximum dynamic range at criticality Q. Such an 
optimum regime typically appears in the bistable regime 
and depends on the past history. 



II. THE MODEL 

As a simple and influential excitable media, we ex- 
plore the Kinouchi-Copelli model 0, [2l[ generalized to 
account for the integration of multiple excitatory inputs. 
We consider N nodes embedded in sparse (Erdos-Renyi) 
random and (Barabasi- Albert) scale-free networks [22| . 
both with an average degree K = 50. Each node i repre- 
sents an excitable unit whose state Si(t) £ {0,1,2} indi- 
cates whether the unit is in the quiescent state [si(t) = 0], 
in the active state [sj(i) = 1], or in the refractory state 
[si(t) = 2], The dynamics obeys probabilistic rules with 
a synchronous update, and St = 1 ms is the discrete time 
step. Every node i at time t updates its state as follows: 



• In the active state Si(t) 
refractory state Si(t + St) 



1, it switches to the 
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• In the refractory state s,(i) = 2, it returns to the 
quiescent state Si(t + St) = with probability p 7 = 



2 




FIG. 1. Continuous and discontinuous spontaneous activity F versus p\. Mean-field approximation (MF) and numerical 
results for random networks of (a) non-integrators (9 = 1), (b) integrators (9 = 2) for both n and Too integration times with 
N = 5, 000; and (c) integrators with different threshold values, for ti and N = 1, 000. Other parameter values are Apx = 0.0025 
and F = 3%. 



• Nodes in the quiescent state Sj(t) = become 
active either (i) by an external driving (or spon- 
taneous activation) with probability ph = 1 — 
exp(— hSt) per time step, where h is the rate of 
a Poisson process; or (ii) by the integration of the 
contributions received from their active neighbors, 
with probability p\ . 

In order to model the integration process, we count the 
number of neighbor contributions Aj(t) received within 
the time window of width t: (t — r,t). In the absence of 
external driving, a node i spikes if Aj reaches at least 9 
inputs, i.e., Aj(t) > 9. Two extreme limits of integration 
time are of particular interest: the infinite integration 
time t — > oo (tqq), where the integration window takes 
into account the entire current quiescent history of the 
node, and a coincidence detection r = 1 ms, (ti), where 
the integration time is limited to St. 



III. CONTINUOUS VERSUS DISCONTINUOUS 
PHASE TRANSITION 

In the absence of external driving (h = 0), the standard 
model without integration (9 = 1) leads to a continuous 
phase transition Q. The average firing rate F, calculated 
over all nodes and over a large time window (10 s), grows 
smoothly for increasing coupling strength above the crit- 
ical value p\ (Fig. Q}i). The critical point is determined 
by the largest eigenvalue of the network adjacency ma- 
trix [l4| ■ For a random network (Fig. QJ,) , the critical 
value is p c x = K , when the average number of spikes 
induced by each spike (branching ratio) is one Q. Con- 
versely, in the presence of integration (9 > 1) the phase 
transition occurs abruptly, generating a bistable phase 
with a hysteresis cycle (see the mean-field approach be- 
low). We calculated the hysteresis cycles by varying p\ 
upward and downward along the whole range in small 
steps of Ap\, activating at each change of p\ a small 
fraction of nodes (Fq, from 1% to 3%) to allow the sys- 
tem to escape from the resting configuration. As shown 




FIG. 2. Dependence of the nature of the phase transition or- 
der on the integration time r for networks (N = 5, 000) com- 
posed of a mixture of both integrators (with a density d of 
9 = 2 nodes) and non-integrators. The solid (dashed) line cor- 
responds to the random (scale-free) network for Ap>> = 0.001 
and Fo = 1%. The left-hand side of the curve corresponds to 
a continuous phase transition whereas the right-hand side cor- 
responds to a discontinuous phase transition. The error bars 
correspond to the standard deviation over ten trials. The 
black open symbol depicts the mean-field shift in the order 
of the phase transition. The left inset panel compares the 
mean-field approximation with the simulations for the den- 
sity of integrators d = 70% and rj. The right inset panel 
illustrates a discontinuous phase transition for r = 2 ms and 
d = 90%. 



in Fig. [TJd, the change in the nature of the phase transi- 
tion is observed for any value of the integration time, as 
well as in the mean-field approximation. The discontinu- 
ous phase transition is also robust for any value of 9 > 1 , 
illustrated in Fig. [TJ: for tj. It can be also seen from the 
figure that larger threshold values generate larger hys- 
teresis cycles. 

While the previous analysis assumes identical nodes, 
next we consider heterogeneous populations composed of 
both integrators (9 = 2) and non- integrators (9 = 1) 
nodes. This situation corresponds to the intermediate 
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FIG. 3. Response curves and dynamic range in networks of integrators (8 = 2) with N = 5, 000. (a) Family of response 
functions for a random network with n and (from right to left) p\ = 0,0.1,0.12,0.14,0.16,0.18. (b) and (d) Dynamic range 
versus coupling strength for (b) random and (d) scale-free networks. At the bistable region, the dotted line (bottom) stands 
for initial conditions with a high activity level and the continuous line (top) stands for initial conditions with a low activity 
level. The error bars correspond to the standard deviation over six trials of two realizations, (c) Dependence of the maximum 
dynamic range A max on the density of integrators in random networks. 



configuration between integrators, as in Fig.[T]D, and non- 
integrators, as in Fig.[]J,. For random and scale-free net- 
works, the minimum density of integrator nodes (d) that 
yields a discontinuous phase transition depends on the 
integration time scale r, as shown in Fig. [2j Although in 
both cases the density of integrators needed to display a 
discontinuous phase transition decreases with increasing 
integration time, the scale-free network requires a lower 
density of integrators. The integration time is funda- 
mental to bind the collective dynamics together. Coin- 
cidence detection restricts the scope of action of the in- 
tegrator nodes and the network is effectively split in two 
parts according to the threshold values. For example, in 
a random network with tj and a density of integrators 
below 80%, the dynamics is dominated by the sub-group 
of active non- integrators, leading to a continuous phase 
transition. In this case of continuous transition, the inte- 
grator nodes do not interfere much in the dynamics: The 
effective connectivity is K(l — d), and the expected criti- 
cal point for the phase transition is given by p\ ~ KI J_ d ^ 
(for the left inset panel of Fig. 2: K = 50, and d = 0.7, 
p c x = Yjr). For larger integration times (r > tj), the dis- 
continuous phase transition (as exemplified by the right 
inset panel) gradually dominates, and the right side of 
the transition increases with r. In this case, the inte- 
grator nodes, although spiking less, tend to remain ac- 
tive, furnishing clear influence in the collective dynamics. 
Therefore, the prevailing dynamics carries the integrators 
finger-print given by the discontinuous phase transition. 



IV. DYNAMIC RANGE 

So far we have analyzed the behavior of the excitable 
media in the absence of external stimuli. In the remain- 
der, we are interested in the response of the system as 
a function of the external driving, considered as a Pois- 
son process with rate h. As illustrated in Fig. the 
response functions for different coupling p\ grow with 
external driving rate h and saturate at a maximum firing 
rate F max = 1 _ 1 ms" 1 , which is determined by the re- 



fractory period p y . Among the response functions, there 
are three regimes. For very low coupling, the response 
functions arc subcritical, the self-sustained solution is 
not allowed, and the activity dies out when h — > 0. On 
the high coupling limit, small perturbations lead the sys- 
tem to the self-sustained mode. In between both regimes 
there is the bistable region. Response functions in this 
regime are history dependent. Very small perturbations 
are typically not enough to drive the system to the self- 
sustained mode. However, at a certain external driving 
rate (which is trial dependent) the system becomes ac- 
tive, as depicted by the upward arrow in Fig. |3Ji. On 
the contrary, if the response function is calculated by re- 
ducing the external driving (leftward arrow), the system 
maintains a high firing rate and the activity does not die 
out when h — > 0. This path dependence could explain 
the large fluctuations in the experimental response func- 
tions as well as the dependence on the measurement time 
period in the olfactory system [23| . 

The bistable regime also confers path dependence to 
the dynamic range. Figure depicts the key elements 
of the standard dynamic range definition. The two hor- 
izontal dashed lines stand for i*b.i (bottom) and fb.g 
(top). They correspond to 10% and 90% of the max- 
imum firing rate (-F max ) subtracted from the minimum 
firing rate [F(h — > 0)], and they cross the response func- 
tions respectively at the external driving intensities of 
ho.i and ho.g. The dynamic range is thus defined as 
the number of decades comprised between ho.i and ho.g- 
A = 10 log Figures [3}d and [3pl show the dynamic 
range for networks of integrators with different integra- 
tion times r, for random and scale-free networks. In the 
bistable regime, when a high firing rate is observed (bot- 
tom line) , the system is only able to distinguish the input 
level intensity. For a low firing rate (top line), the sys- 
tem not only distinguishes the input intensity but also 
detects the abrupt change in the firing rate. The sys- 
tem displays the largest dynamic range in the low firing 
rate and the maximum appears in the bistable regime. 
The height and width of the peak of the dynamic range 
curves depend on the integration time. Coincidence de- 
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tectors show a poor capacity to distinguish the incoming 
input (lower peak and narrower width of A as a function 
of p\). However, for large enough density of integra- 
tors, the dynamic range increases with longer integra- 
tion times (see Fig. , which increases the capacity to 
discriminate incoming inputs. Table U compares the dy- 
namic range of the neuronal networks with and without 
integrators: The maximum enhancement of the dynamic 
range as a consequence of the collective behavior [i.e.. 
A max — A(pa = 0)] is over four times larger in the pres- 
ence (than in the absence) of integration in both random 
and scale-free networks. 





non-integrators (9 = 1) 


integrators (6 


= 2, Too) 


network 


A max _ A ( Q ) 




A max _ A ( Q ) 


y^max 


random 


10 


26 


41 


57 
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7 


23 


32 


48 



TABLE I. Dynamic range (dB) for network size N = 5, 000. 



V. MEAN-FIELD APPROACH 
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FIG. 4. Bifurcation diagram of the mean-field approximation: 
(a) transcritical for 9 = 1; (b) saddle node for 9 = 2. Solid 
lines, stable stationary solutions; dashed lines, unstable ones. 



is given by F t = F t {1) (l -d)+ F^d, where d denotes the 
density of integrator nodes. Then, by generalizing Eq. [T] 
we can find F 1 and F 2 from: 



5 t F^=Q^ Ph + Q^\l-p h )A e t ' 



(2) 



As shown in the left inset panel of Fig. [21 the average 
firing rate of the network qualitatively captures the phase 
transition. 



The mean-field approximation we present corresponds 
to the model version for coincidence detection (tj). For 
K 0, the mean-field map for the average firing rate F 
of a population with threshold 9 can be written as: 

5 t F t+1 = Q tPh + Q t (l - Ph )A e t , (1) 

where Af = [1 — (1 — p\5 t F t ) K ] e is the probability of a 
quiescent node to become active in the next time step 
due to at least 9 neighbor contributions within a single 
time step, Q t = 1 — 5 t F t — Rt is the probability of finding 
a site in the quiescent state, and Rt+i = S t F t + (1 — p-y)Rt 
is the probability of finding a site in the refractory state. 
Iterating the map until convergence we get the solution 
of F in the stationary configuration (t —> oo), which 
is used to compare with the simulations. The numeri- 
cal solutions of Eq. ([1]) for various conditions are shown 
in Figs. [1] and |31 For the population of non-integrators 
(Fig. []Ji) we recover the Kinouchi-Copelli equation 0, 
which describes particularly well the behavior in random 
networks. In the presence of integration, the result cap- 
tures qualitatively the behavior of the phase transitions 
(Fig. [TJd) , the response function, and the dynamic range 
(Fig. A bifurcation analysis reveals some aspects 

of the phase transition as a function of the threshold 9. 
In the absence of input (h = 0), as shown in Fig. 21 
for 9 = 1 there is a transcritical bifurcation; for 9 > 1 a 
saddle-node bifurcation and a stable fixed point at F = 
coexist. 

Analogously, one can also extend the results for het- 
erogeneous populations, as considered in Fig. O At any 
time t, we define, for each subpopulation of threshold 6i, 
Ff* , Rt , and Q\ ^ as the firing rate, and the proba- 
bility of finding a site in the refractory and in the quies- 
cent states, respectively. The firing rate of the network 



VI. SUMMARY AND CONCLUSIONS 

We have studied the collective behavior of an excitable 
media where the units integrate incoming signals [l6j]. 
The presence of a minimum density of integrator nodes 
leads the system to an abrupt phase transition. Discon- 
tinuous transitions have been observed experimentally 
and in threshold models [24j| . in models with adaptive 
interactions [25l [26ll . and in the presence of strong non- 
linear coupling [27] . 

As a consequence of the discontinuous phase transi- 
tion, bistability emerges. In the context of neuroscience, 
bistability is known to play an important role in memory 
maintenance [28|. A bistable regime composed of a con- 
figuration with high or low activity levels [29| has also 
been observed in cortical neurons. Since most neurons 
(if not all) must integrate their incoming post-synaptic 
potentials, our results suggest that the transition to the 
regime of self-sustained activity in a neuronal system 
could be restricted to a discontinuous transition type. 

Concerning the output response to external stimulus 
(which might vary for orders of magnitude) , the bistable 
regime provides two different response types, depending 
on the history (either with low or high activity levels 
for h ~ 0). The low past activity level with an infinite 
integration time gives rise to the largest dynamic range in 
random and scale- free networks. Taking this finding into 
account, biologically inspired artificial stimulus detectors 
with great capabilities can be designed from excitable 
media composed of integrator units [30| . Moreover, we 
expect that our results might also be relevant to other 
systems where integration plays an important role as, 
for instance, in gene regulatory networks fl5j . and social 
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interaction [17[, and it would be interesting to explore 
the behavior of the dynamic range in the recently found 
explosive percolation (3lj . 
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